* ********************************************************************************************
* The Effects of Weather Shocks on Economic Activity: What are the Channels of Impact?
* Sebastian Acevedo, Mico Mrkaic, Natalija Novta, Evgenia Pugacheva, Petia Topalova
* With support from Gavin Asdorian, Olivia Ma, Jilun Xing and Yuan Zeng 
* Replication files for Figure 3
* ********************************************************************************************
clear all
set more off
set matsize 11000

cd "C:\Users\..\replication" // set your working directory to where the replication folder is saved

* the Y variable that will be used in regressions: ngdprpc
global Y ngdprpc

* Horizons
global k 7


* ****************************************************************************************
* Our data - create extra variables
* ****************************************************************************************
use "data/country_dataset.dta", clear
drop if year > 2016

* year squared (for country-specific time trend squared)
gen sqyear = year ^ 2

* create variable for region-year (ry), encode it to be numeric, to be used as fixed effects
gen tempvar1 = wdi_region + " " + string(year) if wdi_region != ""
encode tempvar1, generate(ry)
drop tempvar1

* create country-specific time trend (_IifsXyea*) and time trend squared (_IifsXsqy*)
* dummies for region-year(_Iry*), country (_Iifscode*) and year (_Iyear*)
xi i.ifscode*year i.ifscode*sqyear i.ry i.year


* calculate GDP log difference at different horizons
gen ln_$Y = ln($Y)
by ifscode (year), sort: gen d0_ln_$Y = ln_$Y - l.ln_$Y
forval t = 1/$k {
	by ifscode (year), sort: gen d`t'_ln_$Y = f`t'.ln_$Y - l.ln_$Y
}


* generate squared terms of temperature and of precipitation
foreach X in UDEL_pw10 UDEL_pw90 UDEL_pw50 CRU_pw90 CRU_pw10 CRU_pw50 {
	gen tmp_`X'2 = tmp_`X' ^ 2
	gen precip_`X'2 = precip_`X' ^ 2
}

* ****************************************************************************************
* Variables where to save results
* ****************************************************************************************

gen temperature = _n - 2
gen temperature_coef = .
gen temperature_se = .


* ****************************************************************************************
* Regression (Table 1, column 5, horizon 0)
* ****************************************************************************************
local yy $Y
local tt tmp_CRU_pw50
local pp precip_CRU_pw50
local options "_Iifscode* _Iry*"

reg d0_ln_`yy' `tt' `tt'2 `pp' `pp'2 l.`tt' l.`tt'2 l.`pp' l.`pp'2 l.d0_ln_`yy' `options', cluster(ifscode)
	
* Effect for temperature from -1 to +30 
local row = 1
forval i = -1(1)30 {
	lincom  _b[`tt'] + (2 * _b[`tt'2] * `i')
	replace temperature_coef = r(estimate) in `row'
	replace temperature_se = r(se) in `row'
	local ++row
}	

keep temperature temperature_coef temperature_se
keep in 1/32

tempfile results
save `results', replace


* ********************************************************************************************
* Histograms
* ********************************************************************************************
use "data/country_dataset.dta", clear
keep if year == 2015
drop if tmp_CRU_pw50 == .

gen one = 1

local row 1
gen temperature = .
forval i = -1(1)30 {
	replace temperature = `i' in `row'
	local ++row
}

foreach region in ae em lics {
	gen share_`region' = .
	sum one if `region' == 1
	gen sum_total_`region' = r(N)
	
	local row 1
	forval i = -1(1)30 {
		local j = `i' - 0.5
		local h = `i' + 0.5
		egen sum_t = sum(one) if tmp_CRU_pw50 > `j' & tmp_CRU_pw50 <= `h' & `region' == 1
		replace sum_t = sum_t / sum_total_`region' * 100
		egen sum_t2 = max(sum_t)
		replace share_`region' = sum_t2 in `row'
		replace share_`region' = 0 if share_`region' == .
		drop sum_t sum_t2
		local ++row
	}
}

keep temperature share_ae share_em share_lics
keep in 1/32
merge 1:1 temperature using `results', nogen
order temperature temperature_coef temperature_se

export excel using "output/Figure_3.xlsx", sheet("Stata") sheetreplace firstrow(var)
